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NOMENCLATURE 


Cfj, coefficient for eddy viscosity equation 

f p wall damping function 

k turbulence kinetic energy 

Ip turbulence length scale 

l € turbulence length scale 

Ro rotation number 

Ri Richardson number 

Ri* extra term in equation(4) 

R y turbulent Reynolds number (fc 1//2 y/ v) 

Rt turbulent Reynolds number (k 2 /ue) 

Uj Xj component of mean velocity 
Uj Xj component of velocity fluctuation 
Xj coordinate direction, 1 - streamwise, 2 
e dissipation rate of kinetic energy 
v t turbulent viscosity 

<7fe constant in model equation (3) 

<t € constant in model equation (4) 

a** local mean vorticity 

Cl angular velocity 


transverse 
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ABSTRACT 


The two-layer near wall approach in combination with the k — t model was applied 
to rotating flows. Validation studies show that the conventional length scales have to 
be modified to account for rotating effects. Guided by the Lam and Bremhorst’s low 
Reynolds number model with the rotation-corrected k — e model, a new set of length 
scales formulation is proposed. The developments and validations of the current model 
are reported. A fully elliptic numerical method developed to solve the two-dimensional 
Reynolds- averaged Navier-Stokes equations is employed for the computational study. 
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INTRODUCTION 


The flow fields in typical liquid rocket engines such as the Space Shuttle Main Engine 
(SSME) have the common feature that the flow is highly turbulent and the flow paths 
involve bends through large angles at a small radius of curvature. The flow is normally 
used to drive turbines during its transit through the system. This implies flows with 
significant rotation. Experiments concerning the streamwise curvature effects, as well as 
system rotation effects on turbulence structure have been carried out extensively for the 
last decade[l,2,3,4]. These effects on shear layers are of considerable interests, both as a 
basic fluid mechanics problem as well as a design problem for fluid machineries. 

For these flow problems involving solid boundaries, the turbulence in the low- Reynolds 
number near- wall region has been shown to have dominant effects on the overall flow struc- 
tures. As for the effects of the rotation on the near-wall turbulent structure, Bradshaw[5] 
in 1969 has pointed out the similarities among rotation, streamline curvature, and ther- 
mal stratification on turbulent flows. In all three cases, the extra body force has either a 
stabilizing or destabilizing effect on the flow depending on its interaction with the driving 
mechanism. Near the wall regions, these effects are even more profound. Near the concave 
wall of a curved wall-bounded flow as well as the pressure side wall of a rotating flow, 
turbulence is significantly enhanced with an increase in the production rate of turbulent 
kinetic energy. On the other hand, the turbulence is suppressed near the convex wall and 
even laminarized on the suction side of a rotating flow at high rotational speeds. 

There have been a number of computational studies of predicting the rotation as 
well as streamwise curvature effects on turbulence by using turbulence closure models. 
While some higher-order models including Reynolds Stress models and Algebraic Stress 
Models have been used, majority of the models used are eddy viscosity models based on 
calculations of a characteristic velocity scale and a turbulence characteristics length scale. 



For rotating flows and curved flows, the most common form used is the Monin-Obukhov 
similarity formula (see[ 5] ): 

l 1 

r = ; 777 = 1 - 3 Ri 

lo 1 + l3Ri 

where [3 is a constant, which can be determined from experiments, / 0 is the mixing length 
at zero rotation and in flow with no curvature, and Ri is a Richardson number, defined as 

_ —2 Ro(dU/dy — 2 Ro) 

t= (dUjdy )~ 2 

In which Ro — ClD/U , fl is the angular velocity.The Richardson number is a local stability 
parameter and a negative value of Ri denotes an unstable effect while a positive value 
denotes a stable effect. For wall-bounded flows, these turbulence models have to be used 
in conjunction with models dealing with the near-wall low-Reynolds number region. 

Several procedures have been developed to provide the near wall treatments which 
include using wall functions, coupling to some empirical models, and adding low Reynolds 
correction terms to integrate to the wall. In previous calculations of wall-bounded rotating 
turbulent flows, wall function approaches [6,7] are most frequently used to model the near- 
wall region flow structures. These methods are derived from the logarithmic velocity 
profile based on experimental observations that turbulence at the near-wall region can be 
described in terms of the wall parameters such as wall shearing stress, normal distances, etc. 
Therefore, these methods are not valid if the logarithmic velocity profile no longer prevails 
in the near-wall region. Due to this limited applicability of the wall function methods, 
low Reynolds number models (see, for example [8]) have been developed to overcome 
the shortcomings of the wall functions methods. In Ref.[9], we have extended Chien’s 
[10] low Reynolds model for near- wall turbulence structure calculations, to account for 
rotation/curvature effect. This approach has been shown to be effective in combination 
with added rotation corrections to the k — e model. However, a larger number of grid 
points have to be placed at the near- wall region, and numerical instabilities are frequently 
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encountered. Another alternative is the two-layer model as described in [11,12]. In this 
model, the computational domain is divided into two layers near the wall: the outer layer, 
where the flow is fully turbulent, and the inner layer, where the flow is dominated by 
viscosity and wall damping effects. The turbulent kinetic energy equations are extended 
to include the inner layer, whereas the dissipation rate equations are solved only in the 
outer region. In comparison with the low Reynolds number model, the two-layer model 
usually requires less grid points in the very near-wall region, and generally is more stable. 

The purpose of this study is motivated by the success of the two- layer near wall 
model for separated flows at low speeds [12], and, for high- speed shock/boundary layer 
interactions as reported by Horstman [13]. The approach used here is to apply the two- 
layer near wall treatment to the previously tested rotation-corrected k — e model of [9]. It 
was found that the original length scales formulas used in [11,12] based on non-rotating 
flat plate flow dynamics cannot account for the additional rotation effects in the length 
scale profiles very close to the wall. New formulations are proposed and tested in this 
study. The developments and validations are described in the following sections. 


9 


GOVERNING EQUATIONS 


The Reynolds-averaged momentum equation for incompressible flows on a rotating 


frame can be written as 


dUi 

dxi 


= 0 


dUi d. TTTT . IdP 1 d dUi d-r-r 

dt + dx^ UiU ^ pdxi + Redx^dx^ dxj U ' U] 

(flj Ijfli fljXjflj ) 


( 1 ) 


( 2 ) 


The last two terms in equation (2) are the Coriolis and centrifugal forces, respectively, 
and fii is the angular velocity of the rotating frame. For non-rotating Cases, fl, = 0. 
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TURBULENCE MODELING 


For turbulent flows the Reynolds equations are closed using an eddy viscosity derived 
from the k — t model. The specific turbulence model used here is the one used in our 
previous study involving a correction term to the e - equation. This two equation model 
is expressed as: 


and 


dk d u v.dk. —out 

- + —(Ujk) = —Uv + - « 


de d d ut.dt. 

Ft + dr 1 {U ’ £) -9r, l{, ' + 7; ) d7 i ] 


- C'l -uV — 

1 3 dxj 


-C 2 (l + Ri')j 


Ci = 1.44 


C 2 = 1.92 


(3) 


(4) 


The extra term Ri* that appears in equation (4) is a general Coriolis force modification 
term which has the form 


. k 2 

(c r u;fc — c w n k )n k — 


in which tok = eijk{dUi/dxj) is the local mean vorticity. The model constants C r and 
Co, were fixed to be 0.406 and 8.802, which gave most satisfactory results for all rotation 
numbers tested. 


The two-layer model concept based on Chen and Patel [9] is used here. The computa- 
tional domain is divided into two layers near the wall: the outer layer and the inner layer. 
The turbulent kinetic energy equations are extended to include the inner layer, whereas 
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the dissipation rate equation is solved only in the outer region. Within the inner layer, e 
is specified as 


fc 3 / 2 


e = 


(5) 


The eddy-viscosity within this layer is determined as, 


k 2 


t't 


= C^k 112 !* = f»C» — 


( 6 ) 


where / M = ljl t , the length scale and l ( change exponentially with the turbulence 
Reynolds number, Rk ■ 


Rk 

Ip = Cin[ 1.0 - exp( - )] 


U = Cjn[1.0 - exp( — -^~ )] 


Rk = 


k^ 2 n 


(7) 

( 8 ) 
(9) 


V 


where n is the normal distance from the solid wall. To model the effects of wall junctures 
on the length scale, we follow the procedure commonly used. The resultant length scale 
l e ff is expressed as 

1 m i 

r- = Er (io) 

l 'ff l 3 

where m is the number of walls and lj is the length scale contributed by the j th wall and 
determined by equation (7) and (8). In the present calculations, the values of constants 
are determined as follows: 


Ci = 


At = 


<v /4 

2k 

r 3 /4 


= 2.54 


= 5.08 


( 11 ) 

( 12 ) 


Assigning the value of being 0.95 at the interface of the inner and outer layer, leads 


to 


12 



(13) 


= 70 

The interface of the inner and outer layer is located along a constant grid line where Rk 
is greater than 200. 

Although the two-layer approach has been shown to be effective for representing pres- 
sure gradient effects and provides a more accurate and general description of near-wall 
flow features than computations using wall functions(see[12]), it has not been applied to 
a flow where relaminarization occurs near the wall such as the case for suction side of the 
rotating duct. We first applied the same two-layer model as used for the fully developed 
channel flow to the rotating duct flow and found that it failed to capture the behavior 
of outer regions of the inertia sublayer structure due to the rotation effect. The cause is 
traced back to the formulation of two length scales: i.e. / ^ (which is to be incorporated 
into v t ) and l ( (which is to be incorporated into the sink term of the near- wall fc-equation). 
The original formulations of 1^ and l e were based on flat plate boundary layer formulations 
without rotation effects. To be consistent with the behavior of Coriolis-modified e-equation 
near the wall, we have developed an expression to prescribe l M and l t beneath the matched 
region of the two-layer model. The forms of the l f and are given as 


3U h 0,5 

U = U x [1.0 + 7 (CV— - C w Q)ft(-) 2 ] 

ay e 

(14) 

dU k 1 5 

= u x [ 1.0 + y{c r — - c w D)n(-) 2 ] 

ay e 

(15) 


where 7 = 1.30, These two length scales are then used with Equation ( 5 ) and ( 6 ) to 
compute the eddy viscosity within the viscosity-effected region near the wall. 

Due to the lack of information about length scales near the rotating surfaces, the 
guidance of developing the above formulation has been obtained by the asymptotic match- 
ing analyses as well as the calculated length scale profiles using the low- Reynolds number 
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modeling approach. In the study, the rotation-corrected k — e model ( eq.(3) and (4) ) 
are integrated to the wall with the aid of damping functions of Lam and Bremhorst[14]. 
the advantage of using [14] is that the same model coefficients C i and C 2 are used for the 
two-layer approach. Since the stabilized and destabilized effects are taken care of by the e 
equations, original functions used by [14] are used here. 
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RESULTS AND DISCUSSION 


The governing differential equations were discretized using finite difference method 
based on the control- volume formulation on a non-staggered grid arrangement for all de- 
pendent variables. The velocity-pressure coupling was resolved by the Operator Split- 
ting method (PISO) in time-marching fashion. The detailed methodology and solution 
procedure are described in [15] . In principle, this method is capable of computing two- 
dimensional internal and external flows that are laminar or turbulent, separate or attached, 
incompressible and compressible, steady or unsteady. The method is stable and requires 
no smoothing or explicit under-relaxation other than implied by variations in the time 
step. 

The present model was applied to the fully developing channel flow with and without 
rotations. The computed results will be discussed below. 

Fully Developed Channel Flow 

The calculations imposed a uniform flow profile at the inlet. Turbulent kinetic energy 
was set to 0.003 Uf n and the energy dissipation rate was fixed based on 0.5 channel height 
or step height. A 86 x 61 nonuniform grid was used in calculation domain and about 12 
grid points were allocated inside the near wall layer. The Reynolds number was 33000 
used for the fully developed channel and rotating flows. 

The computational results using the standard k — e model with the two-layer near 
wall model for the velocity, the turbulent kinetic energy, the dissipation rate, the ratio of 
production rate and dissipation rate, the ratio of turbulence viscosity and fluid viscosity, 
the turbulence length scale l ^ and / e , and the wall damping function / M are shown in 
Figure 1 to Figure 7, respectively. All the profiles plotted in these figures were calculated 
using three different grids which are 86 x 47, 86 x 61 and 86 x 81. Among these three grid 
systems, near-wall layer grid points below = 100 is 8, 10 and 12, respectively. It can 
be found that an 86 x 61 grid is sufficient to obtain grid independent solutions. 
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It is noted that the dissipation rate does not vanish in the region very close to the wall. 
The experimentally observed dissipation rate approaches a constant value in this region. 
For y + = 0, the dissipation rate takes the limit value 2uk/y 2 , which is the analytical 
solution of the turbulent kinetic energy equation for a limiting case as y approaches zero. 

In Figure 4, it is seen that the production rate vanishes at the wall and increase to 
the peak value at y + « 15. Hence the assumption of equilibrium turbulent flow may not 
be a good approximation for the region of y Jr = 0 ~ 30. Using the vanishing boundary 
condition for the kinetic energy at the wall yields a growth rate of the turbulent kinetic 
energy and a production rate that are in good agreement with experimental data as well 
as theoretical analysis. For the region of y+ = 30 ~ 100, the equilibrium turbulent flow 
exists since the production rate is approximately equal to the dissipation rate. 

From Figure 5, it can be seen that the ratio of turbulent viscosity and fluid viscosity 
approaches zero at the region y + « 10, which means the flow is laminar. For y + greater 
than 20 the value of the ratio is rapidly increased with linear growth rate. The turbulent 
viscosity is about 35 times fluid viscosity when y + reaches 100. 

The wall damping function shown in Figure 7 is assume to be 0.95 on the interface 
of the inner and outer layer, where Rk is greater than 200. At the wall, the takes the 
limit value of A t j ~ 0.073. 
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Rotating Channel Flow 


For the rotating channel flow, due to different imposed Coriolis accelerations on the 
pressure and suction sides of rotating boundary layers, rotation exhibits destabilization or 
stabilization characteristics on the near wall turbulence. At high rotation numbers, the 
effects on turbulence will lead to eventual laminarization on the suction side and slightly 
enhancement of turbulence intensities on pressure side in rotating channel flow. 

We first present the results using the Lam and Bremhorst low Reynolds near wall 
damping treatments. At solid walls, no slip conditions and k = 0, de/dy = 0 were imposed. 
An 81 x 81 grid was found to be sufficient to obtain grid independent solutions. However, 
it was noted that sufficient number of grids should be packed within the logarithmic layer 
near the wall. 20 grid points were required in our calculations with the first grid points 
away from the wall at about y + =2. 

Figure 8 shows the variation of predicted wall shear velocity with rotation and com- 
parison with the data. The comparison with Johnston’s data shows good agreement with 
the modifying effect on the pressure side of the passage. On the suction side, rotation 
has significant effect on the near wall turbulence depending on the mean flow Reynolds 
number. 

In Figure 9 and 10, the near wall velocity profiles are plotted in the semi- log form 
using the wall parameters. Also on the figures is the logarithmic velocity profile. It is seen 
that fully turbulent flow still exists on the pressure side while laminar flow arises on the 
suction side due to the stabilizing effect for higher rotation numbers. The drastic effects 
of rotation on the near wall turbulence and the law-of- the-wall assumption can be further 
seen from the kinetic energy profiles near both the pressure side wall and suction side wall 
in Figure 11 and 12. The strong destabilization effects of rotation on the suction side of the 
rotating wall cause the reduction of turbulent kinetic energy level. For rotation number 
greater than 0.2, the peaks of turbulent kinetic energy obviously reduce and the turbulent 
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kinetic energy levels are significantly damped. On the other hand, the destabilizing effects 
of rotation on the pressure side slightly enhance the turbulent kinetic energy level. 

To see the rotation effects on the length scale profiles near walls, l € and profiles are 
plotted vs. for Ro = 0.0 and 0.2 at both suction side and pressure side in Figure 13 to 
16. It can be seen that both length scales are drastically effected by rotations. A sequence 
of runs were performed to derive formulations for length scales as a function of rotation 
Richardson numbers described in the modeling section. These formulas are then applied 
for the two-layer approach. 

In the low-layer approach, A 86 x 61 grid was used with 12 grids below the matching 
regions. From the wall to the matching region, the two length scale formulations were 
applied to specify eddy viscosities (Eq.(14)) and dissipation rate distributions(Eq.(15)). 

In Figure 17 and 18 again, the near wall velocity profiles are plotted in the semi-log 
form using the wall parameters. Also on the figures is the logarithmic velocity profile. It 
can been seen that the comparisons are quite good. The effects of rotation on the near 
wall turbulence can again be seen from the kinetic energy profiles near both the pressure 
side wall and suction side wall in Figure 19 and 20. In Figure 21 and 22, the turbulence 
eddy viscosity is plotted for pressure side and suction side, respectively. It shows that the 
effect of rotating on the turbulence eddy viscosity on the pressure side is small, but on 
the suction side the effect significantly decreases the value of eddy viscosity. As rotating 
number approaches 0.2, the eddy viscosity almost becomes zero, which means that only 
molecular viscosity controls the flow field and no turbulent flow exists on the suction side. 

The nondimensional kinetic energy dissipation rate profiles are shown in Figure 23 an 
24 for pressure side and suction side. As rotating number increases there is no major change 
of dissipation rate on the pressure side, but there is a nonmonotonous change existed on 
the suction side. It has been shown that the dissipation rates decrease for small rotating 
number(< 0.1), and increase for large rotating number. According to our experince of 


18 



using low-Reynolds number treatment, such increasing seems to be very sharply and reach 
the much larger when rotating number greater than 0.1. 

The length scale / € and profiles are plotted in Figure 25 to Figure 28. As rotating 
number increases, the length scales do not change much on the pressure side wall. However, 
on the suction side, the length scale l ^ becomes smaller and Z € does not change too much 
as rotation number increases. 
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CONCLUDING REMARKS 

Validation studies of the rotation-corrected k — t turbulence model in combinations 
with a near- wall two-layer model are carried out. A fully developed channel flow was 
first calculated to establish the issues of matching points, matching criteria based on local 
turbulence Reynolds numbers, eddy viscosity ratios and grid point resolutions beneath 
the matching point. Then, the two-layer approach was applied to the rotating channel 
flow case. It was found that the original formula for length scales at the matching point 
of the outer layer and inner layer cannot be applied to the complex near wall turbulence 
structures due to flow system rotation. A new formula for both length scales, guided by the 
Lam and Bremhorst damping function model, was proposed and validated in this study. 
It is shown that with proper corrections in length scale formulations, good predictions of 
near wall turbulence can be obtained. It is recommended that current formula should be 
further generalized for highly curved and separated flows. 
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